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ABSTRACT 



A Steiner Problem in graphs is the problem of finding a set of edges (arcs) with 
minimum total weight which connects a given set of nodes in an edge- weighted graph 
(directed or undirected). This paper develops models for the directed Steiner tree 
problem on graphs. New and old models are examined in terms of their amenability 
to solution schemes based on Lagrangean relaxation. As a result, three algorithms 
are presented and their performance comparred on a number of problems originally 
tested by Beasley (1984, 1987) in the case of undirected graphs. 



RESUME 

Etant donne un graphe G = (V, A) et un sous-ensemble de sommets VI de Y. le 
probleme de Steiner consiste a determiner un graphe partiel de G de longueur mini- 
male permettant de relier entre eux tous les sommets de VI en utilisant eventuellement 
un ou plusieurs sommets de V\\ j. Le present article traite de ce probleme dans le cas 
des graphes non-orientes. Differentes modelisations y sont examinees, ainsi que des 
methodes de solution basees sur la technique de la relaxation Lagrangienne. En guise 
de resultats, trois algorithmes sont presentes et testes sur un ensemble de problemes 
originalement utilises par Beasley( 1984, 1987 ) dans le cas des graphes non-orientes. 



1. INTRODUCTION 



1.1 Graph terminology 

A directed graph G = { \ \ A } consists of a finite set of vertices \ ( | \ = n ) and a 

set of arcs A (\A = m). Each arc in A is an ordered pair of distinct vertices 

v,,Vj (A C U x U). A weighted graph G — {U. A\g}. is the graph G together with 
a nonnegative real function g defined over A. Since our concern is only with weighted 
graphs, we drop the function g from our notation for graphs unless otherwise specified. 

A subgraph G ! C G is a graph G 9 = {U , ,j 4 / }, wdiere l ;l C T and A 1 C A and 
such that each arc a ^ of G 1 ( G A 1 ) is incident to the same vertices in G 1 and in 
G. The weight of a subgraph G 9 , denoted by ir((7 / ), is the sum of the weights c (or 
Cjj for a ^ = (Vj.i'j)) over all the arcs a ^ in A 1 . 

A semipath joining vertices v s and Vj in G is a sequence of vertices v s .v s +\ Vf 

such that the arc (r 5 + ? . j ) £ A or (r s ^ ?4 .i . ) G A for each i. A set of vertices 

U C V is weakly connected if there is a semipath of vertices in U joining any pair 
of vertices in U. A path from v s to Vf in G is a sequence of arcs in A such that 
(tV*WlM r s+l^’ 5 + 2 )i — The path is simple if all the 

vertices {r 5 . r 5 +j, are distinct. A path is a cycle if v s — Vf. A set of vertices 

S C \ r is strongly connected if there is a path from any vertex in 5 to any other vertex 
in 5 which does not contain any vertices in V \ 5. A maximal strongly connected set 
is a strongly connected component of G. An arborcscencc B, (6 C (?) is a set of arcs 
such that 

(i) if (i r u j) are distinct arcs of /?, then ir ? ^ iVj ; (at most one arc 

entering a vertex ) 

(ii) B does not contain a cycle; 

(iii) B is weakly connected. (B is a directed tree.) 

A graph G = {l \ A} is said to be connected if the set \ is weakly connected. 
Let G 9 — {V / ,-4 / } be a subgraph of G = {U .4} and Y" C \ ( % then G* is said to be a 
connected subgraph with respect to \ n if there exists a semipath in G * between each 
pair of distinct vertices in V" . Furthermore, if G 1 is an aborescence, then subgraph 
G * is said to be an aborescence subgraph of G with respect to V r,/ , and if V r/ = V n = 1 
then G 1 is a spanning arborescence of G. Note that one vertex in \ 1 is designated as 
the root for the aborescence. 

1.2 Problem Statement 

Let G = {U, A) be a weighted connected graph. Let \ ) C Y and G(U]) = C 
G: such that G 7 is a connected subgraph with respect to \ \}. The problem is to find 
the least w T eight graph in 6 T (V]). Such graph (denoted by G”(\ i)) is called a directed 
Steiner graph of G with resptd to \\ . If G'(\ \ ) is an aborescence then it is called a 



directed Steiner tree of G with respect to \\ or the minimal arborcscence of G with 
respect to V]. The vertices of the graph 6'*0 i) not from the set V], are referred to 
as Steiner points . 

Our graph notation attempts to tie together the notational conventions developed 
by Bondy and Murty (1976), Hakimi (1971) and Tarjan (1977). An arc from v j to Vj 
is denoted either by (v,*,t»j) or simply by (?, j) if no ambiguity results. 

We restrict the discussion to weighted graphs with strictly positive weight func- 
tions and to the problem of finding the minimal arborescence of G with respect to 
Vj, i.e.. the directed Steiner tree problem. Any solution to this problem, designates 
one vertex in \\ as the root vertex of the resulting minimal arborescence of G with 
respect to 1 \ . To facilitate t he forthcomming mathematical formulations, we aug- 
ment the graph G with an additional vertex denoted as vertex 0 and with a set of 
arcs (O.i’i), for v t 6 The weight function g assigns to these arcs high positive 
identical weights to assure only one such arc in the optimal solution to the directed 
Steiner tree problem. 

The complexity status of the Steiner tree problem on graphs is well settled. Karp, 
1972. proves the NP-completeness of this problem by polynomial transformation from 
Exact Cover by 3-Sets problem. The problem remains NP-complete if all the arc 
weights are equal, if G is a bipartite graph with no arcs joining two vertices in \\ 
(or two vertices in V \ Vj). and also in case G is a planar graph. For more detailed 
reference list on the complexity of this problem see the classic book by Garey and 
Johnson, 1979. As far as a literature review of past research on the Steiner tree 
problem on graphs, the readers are directed to the very fine recent paper by Winter. 
1987. 



2. MATHEMATICAL FORMULATIONS WITH FLOW VARI- 
ABLES 

In the mathematical formulations that follow, the decision variables are of two 
basic types; (i) a binary (selection) variables - j/ z y, which accept value 1 if the arc 
(v t ,Vj) is selected in the solution and the value 0 if not, and (ii) flow variables either 
x j jp (or X(j) representing the amount of flow through the arc (V|,i’y) directed from 
the root vertex 0 to vertex p in \\ (or the total flow on the arc (v n Vj)). 

2.1 Mathematical Formulation - A 

Minimize c,j Viy 

subject to: 



( 1 ) 



— J€t'u{0} x tjp — j€l'u{0} x Jir ~ 



r 1 . for i 0. 

i 0. for / ^ 0./» ^ for all /> G V] ( f >\ 

1 . for i = p, 



VI 


for all (vj.Vj) € -4. p € Vj 


(3) 


o 

Al 


for all (r,, vj) € A.p £ \] 


(4) 


y,j = 0 or ] , 


for all x’j ) 6 A 


(5) 



Given the objective of minimizing the total weight of the selected arcs (the arcs 
for which y t j = 1), the constraints have to ensure the arborescence structure (with 
respect to Vj ). Constraints (2) impose the conservation of flow on all vertices in 

V' \ (Vj U {0}) and ensure that one unit of flow leaves the root vertex 0 for each 

destination vertex p in Ij. Constraints (2) also state that one unit of flow reaches 
each vertex in . Thus, there has to be a path (in G) from the vertex 0 to each one 
of the vertices in . Constraints (3) ensure that flow is allowed only through the 

arcs selected in the solution. Constraints (4) are the flow nonnegativity constraints 

and (5) are the binary value constraints for the appropriate decision variables. 

A nolc: In case the weight function g of the graph G = {V. -4;#} is not strictly 
positive, then in order to ensure a tree structure for the solution to the MDSTP we 
add the constrains = 1 for all j 6 1 j (2a). Since in this paper we restrict 

the discussion to strictly positive functions constraints (2a) are redundent and 
subsequently dropped. 

A more compact mathematical formulation to the one presented in (A) is ob- 
tained by simply aggregating constraints (3) with respect to the index p. In this case 
we obtain: 



x ijp - ^ l y>y for a11 ( v i* v j) € A (3‘) 

We denote the mathematical formulation which contains the equations (1), (2). 
(3'), (4), and (5) as (A 1 ). 

It is of interest to examine the implications of aggregating constraints (3) into the 
form of (3’). Since both formulations construct the same weight solution if solved op- 
timal}', what could be the advantage of attempting to solve one set of equations versus 
the other set of equations ? The answer to this question is in most cases experimental 
in the sense that a number of researchers (Magnanti and Wong, 1981, Cornuejols et 
ah, 1977) were much more successful! in constructing efficient solution schemes for 
the disaggregated formulation (A) than the more compact (A’) formulation. In 
order to formalize this experimental experience we examine the linear programming 
relaxations (constraints (5) become > 0. and y 2 j < 1 for all (v^Vj) G A) of both 



- 5 - 



formulations. Denote the objective function value for the LP relaxation for (A) as 
and for (A’) as -ip- 

Theorem 1 : zf p > z£p. 

Proof : It is clear that any feasible solution for the linear programming relax- 
ation of (A) is also feasible for the linear programming relaxation of (A’). In the 
reverse direction this is not true in general. One can illustrate it by examining a 
equilateral triangle graph with 4 vertices, (one in the center) and 7 directed arcs (see 
Figure 1). The \\ set consists of vertices 1 and 2. 




Figure 1 : Equilateral triangle graph. 

The flow variables 2*031 = 2,2*3]] = 2, -r 3 2 1 — 1*^201 — 1 are clearly feasible in 
the LP relaxation for (A’) and not feasible for the LP relaxation for (A) • 

Given two (different) equivalent integer programming formulations (PI) and (P2) 
for an optimization (minimization) problem. Denote by (RP1) and (RP2) the respec- 
tive formulations obtained by relaxing the integrality constraints. U zpp] > ZRpo 
and not equal in general, then formulation (PI) is called a strong formulation 
with respect to (P2). This concept of strong formulation is stated in Geoffrion, 1979 
and Van Roy, 198ft. 

In our case, we proved that (A) is a strong formulation with respect to (A 1 ).' 

The next mathematical formulation for the minimal arborescence problem with 
respect to \\ requires only two indices for its flow variables. It is based on a similar 
TSP formulation from Gavish and Graves, 1982. 

Mathematical Formulation - B 



Minimize ; = E(,, ; ) € ^ y tJ , 



( 6 ) 
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Subject to: 



—j- 


=i 


x kj - ^”=0 x d 


= -1, 


for all 


A 6 l'i 


(A) 


NT 77 
— J 


-i 


x k j - EIU *ik 


= 0 . 


for all 


A € 1 \ 1 j 


(5) 


x u 


< 


I'llsty, 


for all (v t 


r> v j ) € 


A 


(9) 


x u 


> 


0 . 


for all (r 


i,Vj) € 


A 


(10) 


y*j 


= 


0 or 1 


for all ( v 


€ 


A 


(11) 



The formulations (A) and (B) are equivalent in terms of determining the same 
optimal solution. In terms of the linear programming relaxation, formulation (A) 
is a strong formulation with respect to (B) (the proof of Theorem 1 holds for this 
case too). What is also clear is that any feasible solution (in terms of flow) of the 
LP relaxation for (A’) is also feasible for the LP relaxation of (B) and visa versa. 
Thus, (A’) and (B) are fully equivalent. 

3. MATHEMATICAL FORMULATIONS WITHOUT FLOW VARI- 
ABLES 

3.1 A Note on the Minimal Spanning Arborescence Problem - (MSAP) 

When attempting to solve the MDSTP on a graph, one usually reexamines the 
very similar 'easy 1 problem of constructing a minimal cost spanning arborescence 
(MSAP). In all the formulations of MSAP which are known to the authors, flow 
variables are used (see for example Gavish. 1982) in the same fashion as in the MDSTP 
formulations presented above in Section 2.1. Aneja, 1980, presents a formulation of 
the MDSTP without the flow variables. His is a set covering formulation in which 
the number of constraints grows exponentially with the size of problem instances. 
We present a very simple formulation of the MSAP which does not require flow 
variables nor does the number of constraints grow exponentially with the size of 
problem instances. In this new MSAP formulation we ‘borrow’ a version of TSP 
subtour elimination constraints (Miller. Tucker, and Zemlin, 1960). 

Minimize r = T(, j)e.4 c ijVij 

Subject to 

E? =0 y.->= 1 for all j = 1, 2 , n (12) 

v, - Vj -+■ ny,j < (n — 1) for all i,j el’U {0} 
y tJ - 0 or 1 for all i,j 6 I’U {0} 

where v, and vj are arbitrary real numbers. 

- 7 - 



(13) 

(14) 



Note that those cycle (subtour) elimination constraints (14) are for all i.j 6 
\ U {0} including the ‘artificial* root vertex 0 (contrary to i,j f Q in the TSP 
formulation). Also note that the TSP subtour elimination constraints of the type 
£ IF. — 1 and — 1 f° r ^ V would not be appro- 

priate for this MSAP formulation. 

3.2 Set Covering Formulation - C 

First, we present a set covering type formulation for the MDSTP, modified for 
the directed graph by Wong, 1984, and originally presented by Aneja, 1980, for the 
undirected graph. 

Minimize r = E[iJ)eA c ijVij ( 15 ) 

Subject to: 

H(i j) eA jE x ] y,j > 1 for all A r j C V such that 0 € A’] and A ] fl \\ 0 

(16) 

(A j denotes the complement of A ] in \ ) 

y tJ = 0 or 1 for all (i,j) £ A (17) 

In his solution scheme for formulation C, Aneja. 1980. solves optimaly the linear 
programming relaxation of C. An interesting result due to Wong, 1984. is that 
z^p = z ip(zl z ip — z ip)' linear programming relaxations of the set 

covering formulation C and the flow formulation A. have the same optimal values.) 

The number of constraints in formulation C is exponential in the size of the 
problem. There is too little ‘structure' in this formulation to be usefull in Lagrangean 
relaxation schemes. In order to amend this ‘weakeness* we add a number of redundent 



structural constraints in the next formulation. 

3.3 Modified Set Covering Formulation - D 

Minimize r = Z( t j) € A c ijVt 3 ( 18 ) 

Subject to: 

E” = 0 Vij = 1 for al1 ye V, (19) 

Vi, > Hil - 1 (20) 

u, — Uj -+ ny tJ < (n — 1) for all (i,j) € A (21) 



E(ij)<eO( 7 >) Vij > 1 f° r a11 and all cuts {C'(p)} (22) 

where C(p) is a cut (a subset of arcs in A) between the vertex 0 and the vertex 
P € V'i 



b - 



y,j — 0 or 1 for all ( i.j ) £ A 

v, and Vj are arbitrary real numbers. 



(23) 



The objective function (18) and the constraints (22). (23) define the MDSTP in 
a formulation equivalent to C. Const raints (19), (20). and (21) are redundent. The 
constraints in (19) ensure that only one arc enters each vertex in l j . Constraints (20) 
state that at least C c | — 1 arcs are selected in any solution where initially 1 c = l) and 
in principle \\ C l' c . Constraints (21) are the subtour elimination constraints. Con- 
straints (21) are not the ‘complicating’ constraints in this case. The complicating 
constraints are the set covering constraints (22) which ensure a tree structure solution 
which spans all the vertices in \\. Dropping the constraints (22) results, through the 
(21) constraints in a minimal weight forest solution to the remaining problem. 

Corollary 1: > - C LP {= z^ p ). 

This result simply follows from the fact that we have added constraints (19) to 
an equivalent formulation to C. The other constraints types ((20) and (21)) would 
be satisfied in t lie C formulation through constraints (16). 

4. LAGRANGEAN RELAXATION 

In this section we describe a number of Lagrangean relaxations for the mathe- 
matical formulations presented in Sections 2 and 3. Lagrangean relaxation approach 
for solving ‘hard' problems is based on the observation that by removing the compli- 
cating constraints from a mathematical formulation, the resulting problem is ‘easily* 
solvable. A solution to the relaxed problem constitutes a lower bound on the solution 
to the original problem. The thrust in such an approach is to obtain a maximal lower 
bound which, if it does not solve the original problem, can be integrated into an 
implicit enumeration scheme such as branch and bound. 

4.1 Lagrangean Relaxations of (A) 

We present two Lagrangean relaxations of formulation A. In the first one the arc 
selection constraints are relaxed, resulting in a shortest path type problem. In the 
second relaxation the conservation of flow constraints are moved into the objective 
leading to a more difficult subproblem. 

4.1.1 The First Relaxation 

In formulation (A), the complicating constraints are the arc selection constraints 
(3). which ensure that the unit flow from the root vertex 0 to a vertex p*p € V\ passes 
only through the arcs selected in the solution. Below we present the Lagrangean 
relaxation obtained by moving the constraints (3) multiplied by nonnegative \ lJp 
into the objective function. For a given vector of multipliers A, the problem is 
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Minimize ;(A) — 3I(j t j)gj4 c ijVij — (;j)£.4 ^-p€l’i ^ijp(Vij x ijp ) 

which after rearrangement of terms has the following form 

Minimize c(A) = — ^2p£\\ ^ijp)Vij ~f ^p€'i ^ijp x ijp (24) 

Subject to: 



{ 1, for 1=0, > 

0, for i ^ 0,p > for all p £ Vj 
- 1 , for i = p, j 
0 < ?ijp < 1? for all (i,j) e A,p £ Vj 
y t j = 0 or 1, for all (i.j) € A 



(25) 

(26) 

(27) 



The 4 best* Lagrangean value is obtained by maximizing r(A) over nonnegative 
A's (\ij > 0, for all (i,j) € A.p £ lj). For those optimal Lagrangean multipliers 
(A t jp) the following relationship holds (see Gavish, 1978, and Wong, 1984. for the 
dual formulation of A): 



Cij ~ ZpO\ X ijp ^ 0 for all (i.j) € A (28) 

By observing that the above ((24) - (27)) Lagrangean formulation has the Inte- 
grality Property (Geoffrion, 1974). the lower bound value obtained for the MDSTP by 
solving (24) - (27) is equal to the value for the linear programming relaxation to the 
problem. In this case, the main advantage for examining the Lagrangean relaxation 
of A would depend on how fast, in comparison, can such a relaxation be solved. In 
addition, such a solution could be more amenable for developing good heuristics. 

The (28) inequalities suggest a fast solution procedure to obtain the maximal r(A) 
solulion for the Lagrangean relaxation of A based on the repeated use of the shortest 
path algorithm with cost modifications along the way. In the algorithm outlined 
below, we succssively adjust the values of the A^’s using the subgradient method 
described in Held et al. ( 1974) while preserving the dual feasibility of these multipliers 
via the (28) inequalities. We enforce for all (i,j) € A throughout the algorithm 
below (in each modification of multipliers) the following constraints : Y.p£ V) ^ijp ~ 
c t j. Then the objective of the Lagrangean problem is equivalent to a shortest path 
problem. 

Denote by :( LAl(p)) the value of z(LA \ ) obtained by the Algorithm LAI outlined 
below, given that p, (p £ Fj) is the root node of the Steiner tree. Denote by ~z(LA\(p)) 
the weight of the tree generated by z(LA\(p)) solution (i.e., assign the actual costs 
Cij to the arcs in the Steiner tree). Let Sp{l,p) denote the shortest path from l to p 
in the network with arc costs (A,^). (?, j ) £ .4. 
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Algorithm LAI 

Step 0 : (Initialization) 

Set \ IJp = for all (/,;') € A,p € 1) 

~(^)best ^ ’ her 0. Implter — 0; 

Idivide = (preset parameter); 6 - 2; Iterlim = Maximum number of iterations; 
= smaller gap on b\ei csi = relative precision desired on z(\)i cst \st = 0 

Step 1 : Solve the following (Lagrangean relaxation) problem: 

mine] = Hp£\\ A ijp^ijp subject to constraints (26) and (27). 

This problem can be solved with a shortest path algorithm as follows: 

For each p 6 1) compute all the shortest paths from p to k £ 1) using arc costs 
{A,jp} and denote its cost by Sp. Select the solution with the minimal S p value. 

Step 2 : Updating the bounds 

Lower bound : 

If z{\) > z(X) besi then 

-(A )f, £si = r(A), and for all (i,j) € A.p € 1 j , \ IJp (best ) = A ,jp.~],jp{best) = 
hjp* sl(besi) = si , Implter = 0 

Upper bound : 

Consider the subgraph G 1 — (U U {0},.4C where .4* = {(i*j)l x t jp > 0 for some 
p € V,}. Obtain the shortest path tree structure for this subgraph by computing the 
shortest paths from 0 to each node of Vj. Let {y t] — 1} for each arc (i.j) in this 
shortest path tree and 0 otherwise. 

The current solution vector (y,x) is always feasible for the original problem and 
its cost z f (y,x ) = 'Z( l , J )eA c ijVij 

U Zf(y,*) < -best ihen 

z best = z f(y,r), Implter = 0 
Step 3: Updating A 
(i) If Implter = Idivide then 

Change the value 6 and restart from the best solution z(\)i cst (i.e.. h = h 2, 
Implter = 0, z( A) = z{\) b(st .st = st(best )/2. \ lJp = ^,jp{bcsl),~,, J} , = ‘hjp(bcst) for 
all (/,;') € A.p € lj ). 



ll - 



Else : Compute ihe ascent direction -) and 1 lie step st for the current 

solution. 



1,JP = x ijp ~ yij-. for all (j\j) € A.p € \) 

~-(A)) 

Ihlr 

( i) ) A , jp = A t jp — s t ijp i (i. j) £ A. p £ l ] 

The new multipliers A ,j p are obtained by solving the following problem: 
min{ |A - A||;33 p€ v, A tjr = c tJ , X tJp > 0,(i,j) € A,p £ Vj } 

The vector A is the projection of A on {^pef, A IJp = c,j.A,j p > 0<(i.j) € -4,p € 



Not( : The solution of the projection problem follows the procedure suggested 
by Held et al. (1974. pp. 77). 

Step J, : Sloping Conditions 

(a) If It fr > Iterhm then Stop 

(b) If A < then Stop 

(c) If ~ be \( X ]~' €Sl < ^Ust then Sto P 
Otherwise Go To Step 2 



4.1.2 The Second Relaxation 

Following a different relaxation approach, we dualize on the flow constraints (2) 
using \ ip as the Lagrangean multipliers. The number of multipliers reduces to r?T] . 
which is considerably less than in the previous relaxation. The objective function 
becomes: 

Minimize z( A) - X(,,j) e .4 c ijVij ~ -per, A 0p (^ j€ v x 0 jr ~ Zj£V x jO P ~ U ~ 

—per, A ppi'E.j^v x pjp—'£.jev X jpp+U — — ier,7/o,p— per, A »>(— jer x i jp — — jer x j»p) 

After rearrangement of terms the objective is converted to: 

Minimize r(A) = c ijVij A Hp£\\ — (»,j;)e.4 x ijp^ijp — Uper, ^Opp (30) 

Subject to (3), (4). and (5) where c IJ}) = A Jp — A, p for all (i,j) £ .4 and p 6 1 j. 

Note that the A,j’s in this relaxation are unrestricted in sign and the c ljr can be 
handled implicitly storing n* entries instead of n 3 entries. 
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In order to strenghten the lower bound obtained from such a relaxation we amend 
it with the following constraints: 

y,j > 1 for all j e V] (30) 

and Xjjj = for all (i.j) £ -4 and j G 1] (31) 

The constraints (30), (31) are redundent in the formulation A, but are helpful 
in increasing the Lagrangean bound. Note also that (30) is a relaxation of the (tree) 
constraint which ensures that only one arc enters a node in \j (i.e., \hj — 1 for 

all j £ 1 j ). The ‘tree' constraint is tighter but not easily solvable. 

Before presenting a solution procedure for the problem defined by (29), (3), (4), 
(5), (30), and (31), we make the following observations: 

(i) From the selection constraints (3) we notice that: 

(a) If — 0 then x 1Jp = 0 for all p £ \\ 

(b) If y tJ - 1 then 

(I) x t jj = 1 which follows from (31), and 

(II) x ljp = 1 if c ljp < 0. and x ljp = 0. if c JJp > 0 A/> f j. 

Note that for the same reason as in the first relaxation of A. the maximum lower 
bound for the MDSTP obtained by solving this relaxation can not exceed the bound 
obtained from the linear programming relaxation of the problem. 

Algorithm LA2 

Slip 1: (Initialization) Set X ip = : Iter = 0: lmplter = 0; •(A)j )€5< — ® :z bes1 — 

cost of the best feasible solution found so far: Idivide = preset parameter: t — 2: si = 
0: Iterlim — maximum number of iterations; = smaller gap on — relative 

precision desired on z(\ 

Note that Aj are randomly generated, zi, esi is obtained by computing the shortest 
paths from one node of \\ to all other nodes in Vj. 

Step 2: Solve the Lagrangean Problem: 

Iter = Iter + 1. lmplter = lmplter + 1. Set y tJ — 0, x ljp = 0. for all (i.j) £ 
A , P € l'i. 

( 1 ) For each arc (i.j) € -4 

(i) Compute M l} = c tJ + — rO\,p±j Cjj.) 

(ii) If j € Vj set M,j = M,j -t c tj j 

(iii) If M t j < 0 then set y tJ = 1 and 
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For a]] p € \ ’] set 




.1 for c, Jp < 0 or j = p 



0 otherwise 



(2) For each j € Vj. if M tJ > 0 for all (U j ) £ A then determine / such that 

Mjj = and set yjj = 1, for all p such that c\ jp < 0 or p = j set 

x\ jp = 1, otherwise set x ljp = 0. 

(3) Compute c(A) 

5/cp 3: (Updating the bounds) 

(1) If =( A) > then 

z ( X )best = = *-h P (l>esf) = ^ tp , st(best) = st , Impiter = 0. 

(2) If the solution (y.x) is feasible for the original problem then compute its cost 
and denote ihe value by zp. 

If zp < Zfe si then set zy es1 = zp: Set Impiter = 0; 

Slcp 4 : (Updating the A) 

(i) If Impiter = Idivide then set £ = Impiter = 0. c(A) = z(\)i >(si .st = 
$t(bts1)/2. \,j, = A t p(best = -jtpibest) for all i 6 U.p € Uj. 



f 1st compute the ascent direction */ and the step st 




a - 



In II 



(b ) A jp — A t p st ~iips i £ Up G I] 

Step 5: (Sloping condition) 



(a) If Iter > Iierlimit 



Stop 



(b) If A < 



Stop 




Stop 



< ( best 



- u . 



Otherwise Go To Step 2. 



4.2 Lagrangean Relaxation of (D) 

The complicating constraints in the D formulation are the set covering * (22) 
constraints. The major difficulty is that the number of constraints in (22) is expo- 



tive function multiplied by the appropriate (nonnegative) Lagrangean multipliers we 
obtain the following relaxed formulation: 

Minimize ;(A) = H (l c l} y tJ + Ecyell V r (l - — (i,»6C(p) Vij) 

Subject to constraints (19). (20). (21) and (23) where LI is the index set of all 
cuts Cp between the vertex 0 and the vertices p E \\- 

This formulation can be rewritten as: 

Minimize r(A) = ^ Cr€ n ^ c ~ 

Subject to constraints (19). (20). (21). and (23) where: 



This last problem ((32), (19). (20). (21). and (23)) for fixed values of A can be 
solved in polynomial time. (Again from Lagrangean duality we gel c l} > 0 for all 
(i.j) E A.) For a given A vector denote this problem by LDl. 

The question now reduces to one of finding A' which maximizes the value of ;(A). 
i.e.. :(\* ) = Max x >o{:( A)}. 

Theorem 2 : r(A*) > -f/>- 

Proof : This result is based on the observation that the set of constraints (19), 
(20), (21). and (23) do not have the Integrality Property (Geoffrion, 1974,). 

4.2.1 Dual Ascent Procedure for Initial \ p Values 

In this section we describe a dual ascent procedure for computing a lower bound to 
the undirected version of the MDSTP using the D formulation of the problem. Mod- 
ifing this dual ascent procedure for a general directed graph is left as an algorithmic 
exercise. 



nential in the size of 1 j set. On the other hand the number of y tJ variables is exactly 
n(n • l)/2. This implies that in the linear programming relaxation of this formula- 
tion, most of the constraints in (22) are nonbinding. The difficulty lies in finding the 
binding constraints. By removing the constraints (22) and adding them to the objec- 




€ n A C; . for [i.j) E C(p) for some p E \) 
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First we explain and provide an outline of the algorithm followed by a small 
numerical example and conclude with a detailed description of a dual ascent procedure 
for finding good multipliers to LDl. 

The algorithm begins by solving the problem LDl with — {0} U \\ and all A 
values equal to zero. Let Lq be the cost of this solution. The network structure of this 
solution is that of a ‘sparse' forest F of disconnected components. In case of a directed 
graph, this problem can be solved with a modified Tarjan's algorithm (Tarjan, 1977). 
We differentiate between two component types. Components which contain at least 
one node belonging to \\ U {0} are denoted by T — {7 q, Tj, T 2 , ...}, where 7 q is the 
component which contains the root node. The second set of components consists of 
the points S, 5 C V \ \ \ not contained in any of the components in T. 

At this point we pick a component 7 ^ 0. and compute the minimal cost 

of expanding this component. Let 6 /, = Alin^ a€V\ 77 the multiplier 
value over the cut seperating the nodes in 7^ from all the other nodes. The cost 
matrix is updated to: 



c tJ - 6 k for all j e T k .i 6 V \ T k 
c,j = < (33) 

' c 7J otherwise 

As a result of cost matrix update, one or more arcs on the cut have a reduced 
cost of zero. 

Let Lg ( — Lft t l -f 6^ where t is the iteration number. 

The component 7^ is merged via the zero cost arcs with a number of other 
components. This process continues till 7^ is merged with 7 q at which point another 
component in T is selected. At each iteration at least one component is added to 7). 
thus this merging of components stops at most after \S\ + |T| — 1 steps. When the 
process is completed we have only one component 7 q which contains a subset 5] of 
nodes (Steiner points). S\ C V \ \\. We remove from 7 q all the nodes s,s € S\ which 
are dangling nodes (i.e., their degree is < 1). Let Sq be the set of terminal nodes. 

L B = L B) - 32s€5o c Ps.s 

where p $ is the adjacent node of node s..s € 7o- This last ‘trimming’ step is 
repeated until all the terminal nodes (degree < 1) are nodes in 1 \ only. 

L B is the lower bound value to the MDSTP. 

In order to illustrate those steps we use the following example. The example 
consists of 3 required nodes (1. 3. and 4, where node 1 is designated as the root node) 
and 4 potential Steiner points. The network matrix (the arc weights) is symmetric, 
which in the directed graph version implies two arcs of the same weight and oposite 
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direction between two adjacent nodes of the graph. The initial lower bound obtained 
for this example is 5. The lower bound after two multiplier adjustments is 8. while 
the optimal solution to the problem is 9. 



Example 1 : 
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Figure 2 : The original ‘distance* matrix. 

The initial ‘forest* solution connects node 1 to node 2, node 3 to node 5. and 
node 4 to node 6. This solution has a value of 5. We pick the tree containing the 
nodes 3 and 5 together with the arc from 3 to 5 as our 7* tree in the algorithm. The 
corresponding value is 3 (the arc weight from node 2 to node 3) and the new lower 
bound is 5 3 = 8. The modified distance matrix is as follows: 
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Figure 3 : The modified ‘distance’ matrix after one descent. 

A new tree is selected since the previous one contains the root node 1. The new 
tree is the arc from node 4 to node 6 together with the two nodes. The new 6 ^ 
corresponds to the arc from node 3 to node 4 and has a weight of 2. The new lower 
bound is 5 4* 3 -f 2 — 10. Since the expanded new tree contains the root node and 
there are no trees remaining, we trim the tree from the dangling not required nodes 
and obtain a tree which contains the nodes 1. 2, 3, and 4 and the arcs (1.2), (2.3). 
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and (3.4) for the total of 10-2 = 8 . This solution corresponds to a feasible Steiner 
tree of value of 11 instead of the optimal value 9. 
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Figure 4 : The final modified ‘distance’ matrix. 

4.2.2 The Dual Ascent Algorithm 

Slcp 1 : We start with all A C; = 0 and solve the problem ((32). (19). (20). 
(21), and (23)) to obtain z(0). The solution to this problem is in the form of a not 
necessarily connected set of trees T which might not include all the nodes in Y (i.e.. 
a "sparse* forest). In case we obtain only one tree then we have the optimal solution 
for the Steiner tree problem. Denote by Lg the cost of this solution. 

Step 2 : Let Tj be the number of trees in the forest and let T \ be the set of 
nodes in tree k. One of these trees contains the root node. Denote that tree by Tq. 

IT 

Denote by S the set of nodes in Y \ w ? _q7j (i.e.. the potential new Steiner points). 
Slcp 3 : Pick one of the trees in T \ 7 q. Tree T m for example. Compute 
A = ^in j€Tm ^ Tm {c,j} 

Let: L g = Lg 4- A, c,j = c t j — A. for all i £ T m ,j T m . or j £ T m ,i T m 

Slcp 4 •' ( Merging) Every component in T \ T m and S with exactly one zero 
cost arc to T m (i.e., c t j = 0 j £ T m ,i $ T m ) is merged with T m creating a new 
forest. (In case of components with multiple zero cost arcs to T m see Remark 1.) 
Rename the trees in the new forest. Repeat Steps 3 and 4 until T m is merged with 

To. 



Step 5 : If !T| 7 ^ 1 then go to Step 2, otherwise: eliminate all the terminal 
nodes in Tq which belong to V \ V] and reduce the corresponding Lg value by the 
corresponding arc costs. I.e.. Lg = Lg - Hse 5 0 c p s ,$- Lb is the lower bound value 
for the Steiner tree problem. 

Remark 1 : If we add multiple zero cost arcs between T m and another tree in 
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T then a cycle is created. Thus, at each merge operation only one arc can he added 
between any pair of trees. Since without loss of generality we can assume (different) 
integer arc weights, we can expect that the number of zero cost arcs between a pair 
of trees in a merge step is small. This suggests a parallel processing algorithm for 
the construction of the new trees (merged trees) at Step 5. In each case were more 
than one (say k) zero cost arc exists between a pair of trees, k new merged trees will 
be stored and processed in parallel. At the end, following Step 5, we obtain the best 
lower bound by examining all the trees grown in parallel. 

5. COMPUTATIONAL RESULTS 

The solution methods which were developed and described in 1 he previous sections 
provide a lower and upper bounds on the value of the optimal Steiner tree solution 
for a given graph. Optimal solutions were obtained for the problems for which the 
difference between the value of the upper bound and the lower bound was less than 
one. In order to investigate the comperative performence of the solution methods 
developed in this paper, they were programmed and tested on a set of problems 
taken from Beasley (1984). We present the results of these tests in TABLE I 
below. The algorithms LAI and LA2 were programmed in FORTRAN and the Dual 
Ascent Algorithm (D.A. Algorithm in Table I) was programmed in PASCAL. The 
data set is the one tested in Beasley (1984) and consists of 18 randomly generated 
problems for undirected graphs. For algorithms LAI and LA2 we have considered the 
directed version of these problems by duplicating each arc and assigning directions. 
The three algorithms were tested using a VAX 8600. In all the tests the number of 
subgradiant iterations was restricted to 800. the initial 8 value was set to 2 and the 
parameter Idivide was set to 20 for LAI and to 40 for LA2. 
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Problem 

Number 


in 


1-4 


|ii! 


Algoritl 

Lower 

Bound 


im LAI 
Upper 
Bound 


Algorithm LA2 
Lower Upper 
Bound Bound 


D.A. A1 

Lower 

Bound 


gorithin 

Upper 

Bound 


i 


50 


126 


9 


81.99 


82* 


81.44 


82* 


72 


83 


2 


50 


126 


13 


83.00 


83* 


82.88 


90 


60 


124 


3 


50 


126 


25 


137.88 


138* 


137.07 


177 


107 


160 
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200 
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59.00 


59* 


58.95 
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48 
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200 
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61* 


60.58 


65 


49 
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6 


50 


200 
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121.65 


122* 


121.33 
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188 
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110.97 


111* 


108.79 


123 


92 


138 
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75 


188 
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103.99 


104* 


101.81 


118 
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125 
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75 


188 


38 


219.79 


220* 


207.81 


234 


207 


255 


10 


75 


300 


13 


85.90 


86* 


84.81 


124 


61 


175 


11 


75 


300 


19 


88.00 

| 


88* 


87.79 


127 


64 


192 


12 


75 


300 


38 


172.20 


174 


166.20 


228 


139 


204 


13 


100 


250 


17 


165.00 j 


165* 


162.77 


192 


94 


263 


14 


100 


250 


25 


234.91 


235* 


224.56 


277 


131 


310 


15 


100 


250 


50 


317.60 


318* 


301.27 


353 


249 | 


! 387 ! 


16 


100 


400 


17 


127.00 


127* 


122.65 


162 


73 


256 


17 


100 


400 


25 


128.17 


131 


124.52 


143 


101 


163 


18 


100 


400 


50 


215.56 : 


218 


209.52 ' 


280 


182 


260 



TABLE I : Computational results for ihe three algorithms for Steiner tree 
problem on graphs. 

Out of ihe 18 problems attempted. 15 were solved optimally by the LAI algorithm. 
Only the firsl problem was solved optimally by the LA2 algorithm and the Dual 
Ascent algorithm did not produce a single optimal solution. In terms of the quality 
of the lower bound values, LA 1 s lower bound values dominate the values generated 
by LA2 and the dual ascent algorithm. (The optimal solution is noted by *.) 



SAMMERY AND CONCLUSIONS 

We have presented a number of mathematical formulations for the directed and 
undirected Steiner Iree problem on graphs. These formulations have been used to 
develop Lagrangean based lower bounding procedures for the problem. In compu- 
tational lesls (on 18 problems used by Beasley (1984. 1987) for testing undirected 
Sleiner tree problems), il has been shown that one of the algorithms (LAI) generates 
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lower hound values that are close to the optimal solutions. The nonfeasible solutions 
generated by the Lagrangean based procedure have been incorporated into heuristics 
which attempt to generate “good** feasible solutions. Here again algorithm LAI has 
generated optimal solutions to 15 out of the 18 problems. For the other 3 problems 
the gap between the feasible and lower bound values were under 2 C A. Note that in 
case of an undirected graph. Beasley ( 1984 ) reports solving to optimality only 6 out of 
the 18 problems, and Beasley (1987) again for the undirected graphs reports solving 
to optimality 15 out of the 18 problems using Cray X-Mp/48 machine. The combined 
results (Beasley and ours) solve optimally 17 out of the 18 problems. 

Based on the above results, it is our believe that algorithm LAl can be used as 
an effective tool in B* ranch and Bound based procedures for solving the problem. 



Acknowledgment : We thank J.E. Beasley for providing us with a copy of his test 
problems. 
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